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^ ; Abstract 

First, we have ensured that spherical nonrotating collisionless systems 
\Q . collapse with almost retaining spherical configurations during initial contrac- 

tion phases even if they are allowed to collapse three-dimensionally. Next, on 
^ ' the assumption of spherical symmetry, we examine the evolution of velocity 

dispersions with collapse for the systems which have uniform or power-law 
density profiles with Maxwellian velocity distributions by integrating the 
collisionless Boltzmann equation directly The results show that as far as 
the initial contraction phases are concerned, the radial velocity dispersion 
never grows faster than the tangential velocity dispersion except at small 
radii where the nearly isothermal nature remains, irrespective of the den- 
sity profiles and virial ratios. This implies that velocity anisotropy as an 
initial condition should be a poor indicator for the radial orbit instability. 
The growing behavior of the velocity dispersions is briefly discussed from the 
viewpoint that phase space density is conserved in collisionless systems. 
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1. Introduction 

Since van Albada (1982) demonstrated that collapse of stellar systems 
with small initial virial ratios leads to the end products whose density pro- 
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files are well-described by the de Vaucouleurs _R 1//4 law in projection, much 
attention has been paid to cold dissipationless collapse from the standpoint 
of the formation of elliptical galaxies. In addition, it has been found from col- 
lapse simulations (Polyachenko 1981, 1992; Merritt and Aguilar 1985; Barnes 
et al. 1986; Aguilar and Merrit 1990; Londrillo et al. 1991; Udry 1993) that 
cold enough systems are deformed into triaxial configurations by the radial 
orbit instability (ROI). On the other hand, observations suggest that ellipti- 
cal galaxies are triaxial systems supported by anisotropic velocity dispersions 
(Franx et al. 1991). Putting these things together, we may favor the view 
that the ROI accompanied by cold collapse can naturally explain the ob- 
served properties of elliptical galaxies. 

In spite of the importance of the ROI, its physical mechanism remains 
still uncertain. Many workers have searched for the criterion of the ROI from 
the initial conditions adopted. In particular, velocity anisotropy is considered 
a key ingredient because the ROI could be a kind of Jeans instability in the 
tangential direction as stated by Merritt (1987). Although the ROI appears 
to occur in systems with larger velocity dispersion in the radial direction 
than in the tangential direction, Udry (1993) has reported that for centrally 
concentrated density profiles, triaxial systems formed through the ROI even 
in the initial anisotropic models characterized by 2T ra d/Tj_ <C 1, where T rac j 
and Tj_ are the kinetic energies in radial and tangential motions, respectively. 
In retrospect, velocity anisotropy as an initial condition may have nothing 
to do with the ROI. In fact, numerical studies carried out so far demonstrate 
that collapsing spherical nonrotating systems remain spherical until around 
the maximum collapse times at which the virial ratios become maximum (see 
e.g., Figure 3 of Aguilar and Merritt 1990). Nevertheless, no one has argued 
the velocity anisotropy just before the onset of the ROI. 

These studies mentioned above indicate that we still lack the exact knowl- 
edge of how the velocity dispersions of stellar systems, initially far from equi- 
librium, evolve with collapse. Recently, Kan-ya et al. (1995) have shown that 
for spherically symmetric stellar systems the tangential velocity dispersion 
grows faster than the radial one as far as the initial contraction phases are 
concerned. Unfortunately, however, their analyses are based on such a cold 
approximation that the system has nearly zero velocity dispersions initially. 
In reality, stellar systems would form with some degree of velocity dispersion. 
If the system with a power-law density profile has velocity dispersion at its 
birth, individual mass shells can overlap with one another during the collapse, 
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which is then likely to affect the growing behavior of the velocity dispersions. 
Aside from the cold approximation, they analyzed the systems on the basis 
of spherical symmetry. Consequently, we cannot apply their results directly 
to real galaxies which would have collapsed three-dimensionally. Therefore, 
we need to unravel the collapsing behavior of stellar systems at early stages, 
and then investigate the growing behavior of the velocity dispersions which 
they have from the beginning. 

First, in section 2, we ensure that spherical stellar systems evolve without 
deviating from spherical configurations at the early stages of gravitational 
collapse. Next, in section 3, on the assumption of spherical symmetry, we 
study the growing behavior of velocity dispersions for the systems which 
have uniform or power-law density profiles. Discussion and conclusions are 
presented in section 4. 

2. Collapsing Behavior of Spherical Systems 

We first study how spherical stellar systems, initially far from equilibrium, 
collapse three-dimensionally at the early stages of contraction. Above all, the 
time evolution of axis ratios, b/a and c/a, is presented, where a, b, and c are 
the major, intermediate, and minor axes, respectively. 

We use the self-consistent field (SCF) method, termed by Hernquist and 
Ostriker (1992), to study the collapsing behavior of spherical stellar systems. 
Since the SCF method requires no introduction of a softening length, it is 
suitable to cold collapse simulations in which the modification of a central 
force field through a softening length can have a large effect on the global 
dynamics. In addition, the reliability of the SCF method is demonstrated 
by Hozumi and Hernquist (1995) who applied it to spherically symmetric 
systems. In the present simulations, we no longer restrict ourselves to spher- 
ical symmetry. Thus, the angular dependence of the density and potential 
as well as the radial one is expanded in the SCF code. Since collapse sim- 
ulations reveal that density profiles of the end products resemble galaxies 
obeying the de Vaucouleurs i? 1 / 4 law in projection, we choose the basis set 
proposed by Hernquist and Ostriker (1992) on the ground that its lowest 
order members accurately describe spheroidal objects like elliptical galax- 
ies (Hernquist 1990). We adopt n max = 16 and / max = m max = 6, where 
n ma x is the maximum number of radial expansion coefficients and Z max and 
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m max are the maximum numbers of angular expansion coefficients. We em- 
ploy N = f 00, 000 particles of equal mass. The equations of motion are 
integrated in Cartesian coordinates using a time-centered leapfrog algorithm 
(e.g., Press et al. 1986). 

The initial conditions consist of those spherical systems of total mass M 
and radius Rq which have density profiles of p(r) oc r~ n , where n is chosen to 
be and 1. The velocity distributions are Maxwellian with the dispersions 
determined from the virial ratio, |2T/W|=10 -:L5 , where T and W denote the 
kinetic and gravitational energies, respectively, although the velocities are 
assigned so that the initial data do not include escapers. The unit of mass 
and the gravitational constant, G, are taken so that M — 1 and G — 1. The 
unit of length is determined from the relation such that \2T/W\xRq = 1. 
Then, we choose the length scale of the basis functions to be 2. 
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Fig. 1. Evolution of the axis ratios for (a) the uniform density sphere and 
(b) the p oc r _1 model. The virial ratios are chosen to be 10~ 15 . The solid 
lines show b/a and the dashed lines c/a, where a, b, and c denote the major, 
intermediate, and minor axes, respectively. The vertical arrows indicate the 
maximum collapse times. 

In figures la and lb, we present the time evolution of the axis ratios, b/a 
and c/a for n = and 1, respectively. The axis ratios are computed from the 
principal moments of inertia for the most bound 80% of the particles. We stop 
the simulations just after the maximum collapse times because we pay atten- 
tion only to the early stages of gravitational collapse. According to Cannizzo 
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and Hollister (1992), the final elongation becomes larger as the initial density 
profile is more centrally concentrated provided that the initial virial ratio is 
fixed. However, as far as the initial contraction phases are concerned, we 
can see that the systems approach the maximum collapse times with almost 
retaining spherical configurations, irrespective of the initial density profiles. 
In fact, figure lb shows that the axis ratios change little from the initial val- 
ues of b/a = 0.997 and c/a = 0.993 during the early contraction times (from 
t = to t ks 170) even for the n — 1 density profile model. The small devia- 
tions of b/a and c/a from unity before the maximum collapse phases are, in 
large measure, considered root-mean-square fluctuations coming from finite 
numbers of particles used in the simulations. We conclude, therefore, that 
spherical stellar systems which have uniform or power-law density profiles 
do not deviate considerably from spherical configurations during the early 
stages of contraction even if they are allowed to collapse three-dimensionally. 

3. Growing Behavior of Velocity Dispersions 

We have demonstrated in section 2 that spherical stellar systems collapse 
symmetrically in configuration to a considerable degree at the early stages of 
contraction. Thus, spherical symmetry can be a good assumption when we 
study the growing behavior of velocity dispersions for such systems as long 
as their evolution stays in the initial contraction phases. 

3.1. Analytical Predictions 

Before proceeding into numerical simulations, we summarize the growing 
behavior of the velocity dispersions predicted from analytical results. In 
figure 2, we reproduce the evolution of the velocity dispersions for spherical 
symmetric systems on the basis of Kan-ya et al.'s (1995) analyses. It can be 
noticed from figure 2a that for the uniform density sphere the radial velocity 
dispersion grows exactly in the same manner as the tangential one and that 
the growth of both the dispersions is independent of radius. Figures 2b and 
2c show that for the power-law density spheres the radial velocity dispersions 
decrease after some collapse time, while the tangential velocity dispersions 
continue to increase with time at all radii. In particular, the radial velocity 
dispersions decrease inward radially after some collapse time for the p oc r~ 0m5 
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model, and always decrease inward radially for the p oc r _1 model. Of 
course, this behavior of the radial velocity dispersion comes from Kan-ya et 
al.'s (1995) assumption that the system can continue to contract infinitely. 
Furthermore, they neglect the existence of the outer edges. Therefore, their 
analyses would apply to some intermediate radii of the systems which have 
finite extent and some degree of velocity dispersion from the beginning. 
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Fig. 2. Evolution of the velocity dispersions on the assumption that the 
effects of velocity dispersions are small enough, for (a) the uniform density 
sphere, (b) the p oc r~ 0,5 model, and (c) the p oc r _1 model. The solid lines 
show the square of the radial velocity dispersion of, and the dashed lines 
that of the tangential velocity dispersion of . Note that the evolution of the 
radial velocity dispersion is exactly the same as that of the tangential velocity 
dispersion for the uniform density sphere. The velocity dispersions and the 
radius are normalized with the initial value <7q and the cutoff radius Rq, 
respectively. The circles are added on the o\ lines to indicate the evolutionary 
stages of the collapse: the filled circle is the half mass radius and the open 
circle the radius within which the system contains 10 % of the total mass. 
The time t c is the collapse time at the cutoff radius Ro. 

In summary, we can see that the growing rate of the radial velocity dis- 
persion never exceeds that of the tangential one with a cold approximation. 
In particular, for the power-law density spheres, the ratio of the radial to 
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tangential velocity dispersions becomes smaller as time proceeds. In addi- 
tion, the decrease of the ratio is accelerated as the density gradient becomes 
steeper. 

3.2. Numerical Results 

Now that spherical symmetry is adopted, we use a phase-space method de- 
veloped by Fujiwara (1983) because it gives desirable results insusceptible to 
random fluctuations. This method enables us to obtain a smooth distribution 
of velocity dispersions along radius. We employ (N r , N u , Nj) = (150, 151, 80), 
where N r ,N u , and Nj are the numbers of mesh points along radius, radial 
velocity, and angular momentum, respectively. The units of mass, M, and 
the gravitational constant, G, are taken so that M = 1 and G — 1. The unit 
of length is determined from the relation such that \2T/W\xR = 1 as in 
section 2. 

We examine the time evolution of the velocity dispersions by gravitational 
collapse for density profiles of p(r) oc r~ n with Maxwellian velocity distri- 
butions, where n=0, 1/2, and 1 are chosen, on the assumption of spherical 
symmetry. The initial virial ratios are taken to be 10 -1 ' 5 , 10 _1 , and 10 -0 ' 5 
for each value of n. 

In figure 3, the time evolution of of and of is shown for the uniform 
density spheres, i.e., n=0, where o r and o t are the radial and tangential 
velocity dispersions, respectively. It can be noticed that the evolution of o^ 
and al is very similar to one another, regardless of the initial virial ratios. 
That is, both the velocity dispersions grow at the same rate within some 
radius, which means that the systems keep the isothermal nature there as 
indicated by the analytical predictions. The outer regions, however, show 
the deviation from them. This is because the velocity dispersions which the 
system has from the beginning operate to blur the outer edge: the spread 
of the radial velocity becomes smaller outward to the edge. In this way, the 
isothermal nature disappears with time from outside. This effect is more 
prominent for the larger initial virial ratios. Therefore, the coldest model 
shows the abrupt drop in the radial velocity dispersion (figure 3a), while the 
other models show the gentle decrease of it (figures 3b and 3c). 
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Fig. 3. Evolution of the velocity dispersions for the uniform density spheres 
with initial virial ratios (a) \2T/W\= 10" 15 , (b) \2T/W\= 10" 1 , and (c) 
\2T/W\= 10~ 0,5 . The solid lines show the square of the radial velocity 
dispersion of and the dashed lines that of the tangential velocity dispersion 
of. The velocity dispersions and the radius are normalized with the initial 
value of and the cutoff radius Ro, respectively. The filled circle is the half 
mass radius and the open circle the radius within which the system contains 
10 % of the total mass. The time t c is the collapse time at the cutoff radius 
Rq with no initial velocity dispersion. 

In figures 4 and 5, we present the time evolution of of and a\ for the 
p oc r~ - 5 and p oc r _1 models, respectively. Now, let us focus on figures 4a 
and 5a. If we pay attention to the intermediate radii, the tendency of the evo- 
lution of both the velocity dispersions agrees roughly with that predicted by 
the analytical results (see figures 2b and 2c). In the central region, however, 
the nearly isothermal evolution can be seen, which differs from the analyti- 
cal results. This is because the velocity dispersions prevent the system from 
contracting infinitely. In fact, the density near the central region was kept 
approximately uniform throughout the contracting phase. Thus, the evolu- 
tion of the central part resembles that of the uniform density spheres: both 
the radial and tangential velocity dispersions increase in the same manner 
as the system collapses. The deviation from the analytical predictions also 
appears near the outer edges. The same reason as described for the uniform 
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Fig. 4. Evolution of the velocity dispersions for the par 5 models with 
initial virial ratios (a) \2T/W\= 1CT 1 - 5 , (b) \2T/W\= 10" 1 , and (c) \2T/W\= 
XQ-o.5 no t a tions and their meanings are the same as in figure 2. 




Fig. 5. Evolution of the velocity dispersions for the par 1 models with 
initial virial ratios (a) \2T/W\= lO" 1 ' 5 , (b) \2T/W\= 10" 1 , and (c) \2T/W\= 
1Q-0.5 ^p ne no t a tions and their meanings are the same as in figure 2. 
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density spheres can be applied to the decrease of the radial velocity dispersion 
near the edges. The larger amount of velocity dispersion does not change 
largely the evolution of of and of as can be seen from figures 4b, 4c, 5b, and 
5c. 

4. Discussion and Conclusions 

We have found that spherical nonrotating systems keep their spherical 
configurations during the collapse at the early stages of contraction. We thus 
carried out collapse simulations on the assumption of spherical symmetry to 
examine how the velocity dispersions of spherical stellar systems evolve with 
collapse. Then, we have demonstrated that the radial velocity dispersion 
never grows faster than the tangential velocity dispersion except at small radii 
where the nearly isothermal nature is retained up to around the maximum 
collapse times, irrespective of the initial density profiles and virial ratios. 

Our results are basically consistent with those of Kan-ya et al. (1995), 
although they assume that initial systems have nearly zero velocity disper- 
sions. Thus, we can see that in general the system just before the onset of 
the ROI has larger velocity dispersion in the tangential direction than in the 
radial direction, which might be different from the situation that we envis- 
age for the occurrence of the ROI. However, these findings can be understood 
qualitatively with relative ease as follows. If the initial density distribution 
is uniform, the collapse time of the system is independent of radius, so that 
the conservation of phase space density results in or x R = const., where 
or and R denote the characteristic radial velocity dispersion and radius of 
the system, respectively. In addition, the conservation of angular momen- 
tum of each star leads to o T x R = const., where ot is the characteristic 
one-dimensional tangential velocity dispersion. Then, the ratio, or/ot, re- 
mains unity during the early stages of contraction when the system starts 
with an isotropic velocity distribution. On the other hand, if the initial den- 
sity distribution is power-law such as p oc r~ n (n > 0), the stars of inner 
parts fall first and then those of outer parts follow. Consequently, the width 
of the radial velocity spread for the power-law density profiles is narrower 
than that for the uniform density profiles while ot x R = const, still holds. 
It may be helpful to remind that the collapse time of the outer edge for the 
power-law density spheres is the same as that for the uniform density spheres. 
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As a result, the ratio, or/ot, continues to decrease from its initial value of 
unity with time if the system has a power-law density profile and an isotropic 
velocity distribution at the beginning. 

In the numerical simulations, the deviation of the evolution of the velocity 
dispersions from the analytical results emerges at small radii for the power- 
law density profiles and near the edge of the system for the uniform and 
power-law density profiles. Since the system has finite extent, the velocity 
difference of the phase space elements makes the radial velocity dispersion 
smaller near the edge as collapse proceeds. On the other hand, the veloc- 
ity dispersions can save the system from contracting infinitely. Then, the 
conservation of phase space density again means the increase of the radial 
velocity dispersion as the system contracts. Concerning the isothermal na- 
ture at small radii, we point out that the density near the central region was 
kept approximately uniform with collapse. 

We can now explain why the collapse simulations thus far have shown that 
spherical systems approached the maximum collapse phases with retaining 
their initial shapes. This is because the tangential velocity dispersion remains 
equal to or overwhelms the radial one as collapse proceeds. Conversely, the 
predominance of the growth of the tangential velocity dispersion prevents the 
system from being deformed into elongated configurations at early collapse 
phases. It is thus suggested that the onset of the ROI is determined from 
the physical conditions at least just before or around the maximum collapse 
phase. Even though we start with the condition that 2T rad /T± > 1, the 
growth of the tangential velocity dispersion can easily catch up with or exceed 
that of the radial velocity dispersion if the system has a power-law density 
distribution. Since we have found that the criterion of the ROI will not 
depend strongly on the initial velocity anisotropy, it is no wonder that small 
values such that 2T Tad /T± <C 1 lead to triaxial end products as reported by 
Udry (1993). Therefore, we conclude that initial anisotropy parameters like 
2T rad /T L should be a poor indicator for the ROI. 

We are grateful to Professors S. Kato and S. Inagaki for their careful 
reading of the manuscript and helpful comments on its contents. 
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